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SUMMARY 

The  results  of  a  new  multichannel  deconvolution  method,  applied  to  array  recordings  and 
three  component  station  networks  for  teleseismic  P  waves,  are  presented  and  interpreted  in 
terms  of  possible  surface  reflections  and  other  arrivals  from  the  Nevada  Test  Site  (NTS), 
Novaya  Zemlya,  and  Eastern  Kazakh  Test  Site  (EKTS)  nuclear  explosions.  The  deconvolution 
method  utilizes  the  well  known  fact  that  P  wave  spectra  can  be  decomposed  into  source  and 
receiver  spectral  factors.  The  source  functions  obtained  in  the  deconvolution  process  provide  a 
better  picture  of  the  nature  of  explosion  source  time  functions,  and  in  particular  of  the  presence 
or  lack  of  secondary  arrivals  following  the  P  wave  such  as  pP  or  spall.  The  presence  of  such 
secondary  arrivals  and  their  effects  on  the  first  cycle  of  the  P  wave  are  very  important  in  yield 
estimation.  For  most  events  at  the  eastern  part  of  EKTS  the  source  time  functions  appear  to 
contain  a  pP  arrival  but  they  also  show  later,  unexplained  arrivals  and  other  complexities.  At 
other  test  sites  often  there  are  no  clearly  identifiable  pP  phases  in  the  deconvolved  traces.  Joint 
deconvolution  of  central  EKTS  data  using  all  AWRE  arrays  indicated  strong  azimuthal  asym¬ 
metries  in  the  body  wave  radiation.  Deconvolutions  of  NTS  events  were  considerably 
degraded  by  the  limited  signal  bandwidth  due  to  strong  mantle  attenuation  under  this  test  site. 
The  site  functions  are  also  complex  in  most  cases.  Site  and  source  effects  contribute  about 
equally  to  the  energy  observed  in  the  P  codas  of  the  events  analyzed.  During  the  second  year 
of  this  project,  we  will  deconvolve  several  additional  data  sets,  apply  the  deconvolution  tech¬ 
nique  to  three-component  data,  use  the  amplitude  and  area  of  the  deconvolved  source  time 
functions  to  develop  improved  estimates  of  yield,  and  assess  the  combined  influence  of  Q  and 
pP  on  yield  at  the  various  test  sites. 
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INTRODUCTION 


The  determination  of  the  effect  of  surface  reflection  and  other  secondary  arrivals  on  the 
estimated  yields  of  nuclear  explosions  has  been  a  major  problem  in  nuclear  monitoring  work. 
The  major  philosophical  issues  in  such  work  are  the  following: 


a)  Are  the  explosion  source  functions  simple  superpositions  of  P  and  pP  pulses  or  are 
more  complex  source  time  histories  appropriate?  Much  of  the  previous  work  on  P 
waveforms  simply  assumed  the  validity  of  the  P+pP  models  and  proceeded  to  fit  the 
model  to  determine  pP  amplitudes  and  delay  times. 

b)  Can  we  predict  the  effect  of  the  secondary  arrivals  on  the  yields  estimated  from  time 
domain  measurements  on  the  first  cycle  of  the  P  wave  signal?  This  depends  on  the 
resolvability  of  secondary  arrivals  and  their  physical  interpretation. 

c)  What  is  the  effect  of  these  secondary  arrivals  on  frequency  domain  measures  of 
yield? 

d)  What  are  the  limitations  imposed  by  Q  and  the  noise  on  the  time  resolution  in 
separating  multiple  arrivals?  Clearly,  it  is  not  possible  to  resolve  pP  delay  times  of  a 
few  tenths  of  a  second  along  low  Q  paths. 

e)  And,  in  general,  what  processes,  e.g.  Q,  scattering,  near  source  and  receiver  wave 
conversions,  spall,  Rayleigh-to-P  conversions,  etc.,  determine  the  waveforms  and 
spectra  of  teleseismic  P  waves? 


In  this  progress  report  we  present  results  of  the  application  of  a  new  maximum-likelihood 
iterative  deconvolution  technique  which  is  designed  to  help  resolve  the  above  questions.  The 
method  is  based  on  the  factoring  of  P  wave  spectral  matrices,  consisting  of  recordings  of  mul¬ 
tiple  events  at  common  test  sites  and  using  a  common  set  of  array  sensors,  into  site  and  source 
related  factors.  The  method  has  its  roots  in  several  previously  proposed  deconvolution  tech¬ 
niques. 

First,  we  discuss  previous  work  in  estimating  pP  times  and  amplitudes  and  previously 
developed  deconvolution  methods.  Then  we  outline  the  maximum  likelihood  multichannel 
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deconvolution  algorithm  and  present  the  results  of  deconvolutions  of  data  from  nuclear  explo¬ 
sions  at  Pahute  Mesa  and  Yucca  Flats  at  the  Nevada  Test  Site  (NTS),  at  the  Soviet  Eastern 
Kazakh  Test  Site  (EKTS),  Novaya  Zemlya,  and  the  French  Sahara  test  sites.  The  deconvolu¬ 
tions  are  interpreted  in  terms  of  the  presence  and  characteristics  of  pP,  the  presence  and  com¬ 
plexity  of  other  secondary  arrivals  and  of  the  P  wave  coda,  and  azimuthal  variations  in  the 
source-time  function.  We  find  that  explosion  source  time  functions,  while  they  generally 
appear  to  contain  a  phase  that  is  readily  ascribable  to  pP,  also  contain  a  number  of  other 
arrivals.  Both  source  and  site  terms  appear  to  contribute  to  the  presence  of  P  coda. 
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BACKGROUND 

One  of  the  fundamental  problems  that  motivate  research  on  pP  is  the  possible  effects  of 
such  secondary  arrivals  on  estimated  yields.  Assuming  that  mb  is  used  for  deriving  yields 
through  some  mb-yield  formula,  a  complex  interaction  of  Q,  source  scaling,  and  pP  interference 
will  affect  the  measured  mb  which  can  be  modeled  by  synthetic  simulation  of  the  waveforms. 
For  this  we  need  reliable  estimates  of  Q,  source  waveforms,  and  pP  parameters  derived  from 
the  data. 

In  much  of  the  early  work  it  was  assumed  that  P  waves  are  composed  of  P  and  pP  pulses 
possibly  convolved  with  some  site-related  transfer  function.  Therefore,  a  variety  of  methods 
have  been  proposed  for  estimating  pP  times  and  amplitudes.  Spectral  methods  (Cohen  1970, 
1975;  Shumway  and  Blandford  1978)  attempt  to  make  use  of  the  spectral  modulation  due  to 
pP.  In  the  case  of  complex  site  related  modulation  of  P  wave  spectra,  the  spectral  minima  due 
to  any,  even  solitary,  pP  arrival  besides  the  direct  P,  may  be  completely  disguised,  ’hus  making 
the  interpretation  of  spectral  minima  difficult.  Interpretation  of  power  spectra  in  terms  of  mul¬ 
tiple  arrivals  is  totally  nonunique  since  the  phase  information  is  lost  and  an  infinite  number  of 
pulse  sequences,  derivable  from  each  other  by  dispersive,  all-pass  filters  (Robinson  and  Treitel 
1980),  may  yield  identical  power  spectra. 

Direct  modeling  of  narrow  frequency  band  P  waveforms  was  also  frequently  used  for  the 
determination  of  various  source  and  path  parameters  in  seismology  (Burdick  and  Helmberger 
1979,  Bache  1982).  It  appears,  however,  that  too  high  values  of  the  attenuation  parameter  t*P 
were  used  in  much  of  the  past  work  and  the  results  are  not  acceptable.  This  was  also  pointed 
out  by  Douglas  (1986).  Moreover,  it  has  been  shown  by  Cormier  (1982)  that  in  the  case  of 
explosion  P  waves,  such  waveform  fits  are  highly  nonunique  because  the  attenuation  parameter 
t*,  pP  times  and  amplitudes,  and  the  various  source  related  parameters  strongly  trade  off 
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against  each  other.  Generally  speaking,  waveform  fits  are  not  satisfactory  substitutes  for 
deconvolutions  which  make  use  of  the  total  information  in  seismic  signals  more  efficiently  and 
take  advantage  of  the  high  dynamic  range  and  consequendy,  the  broad  signal  bandwidth  pro¬ 
vided  by  modern  digital  recording  systems. 

Another  approach  that  has  been  proposed  involves  cross-equalization  (or  as  the  authors 
call  it,  "intercorrelation")  of  narrow  band  short  period  waveforms  (Mellman  and  Kaufman 
1981;  Lay  et  al  1984  a,b).  This  method  produces  "pseudo-pP"  delay  times  and  amplitudes,  i.e. 
the  best  fits  to  simple  P+pP  models.  Since  the  results  are  not  deconvolutions  in  the  true  sense 
they  do  not  prove  the  existence  or  correctness  of  such  models.  (Lay  et  al  1984  a,b;  Mellman 
et  al  1984).  However,  it  is  easy  to  prove  by  simple  simulations  that  very  good  cross¬ 
equalization  (intercorrelation)  of  pairs  of  P  wavetrains  is  possible  with  narrow  band  data,  even 
if  the  wavetrains  are  considerably  more  complex  than  P+pP  combinations.  Moreover, 
interevent  spectral  amplitude  ratios,  which  should  be  of  the  form 


.  -icon, 

1  -  otie 

,  -im, 

1  -  a2e 


(l) 


(because  the  site  factors  should  cancel),  cannot  be  fitted  by  any  reasonable  values  of  aj  <X2,  Xj 
or  x2,  where  the  a’s  and  x’s  are  the  pP  amplitudes  and  delay  times  respectively.  The  forms  of 
the  spectral  ratios  computed  from  data  are  quite  complex,  indicating  complex  sequences  of 
arrivals  from  the  source  regions  (Der  et  al  1986b). 

As  a  general  comment,  it  may  be  stated  that  all  proposed  techniques  seem  to  work  quite 
well  with  synthetic  data  that  conform  to  the  proposed  signal  model.  The  problem  is  that  it  is 
rarely  attempted  to  verify  from  the  data  that  the  underlying  models  are  valid.  Based  on  present 
knowledge,  however,  we  cannot  even  be  sure  that  pP  arrivals  exist  in  most  P  waves  from 
nuclear  explosions.  For  example,  the  spectral  modulations  introducing  pP  "nulls"  in  the  spectra 
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are  very  rarely  observed  and  often  the  associated  "pF’  lag  times  are  unreasonably  long.  More¬ 
over,  it  seems  certain,  that  along  some  low  Q  paths  a  short  pP  delay  time  is  not  resolvable,  and 
thus  no  usable  pP  "correction"  may  be  derivable  from  the  data. 

Time  domain  methods  for  estimating  pP  times  and  amplitudes  include  various  deconvolu¬ 
tion  schemes  (Marshall  1972;  Bakun  and  Johnson  1973;  Kemerait  and  Sutton  1982;  Israelson 
1983)  which  generally  preserve  the  phase  information  and  do  not  depend  on  die  simple  P+pP 
model.  Least  squares  (12  norm)  deconvolution  models  work  well,  but  those  that  use  the  lj 
norm  (Levy  and  Fullagar  1981;  Mellman  et  al  1984)  were  shown  by  Jurkevics  and  Wiggins 
(1984)  to  result  in  erroneous  deconvolutions  as  judged  by  synthetic  simulations.  Deconvolu¬ 
tion  methods  belong  to  the  general  type  of  "appraisal"  methods,  borrowing  the  terminology  of 
Oldenburg  (1981),  in  which  we  attempt  to  establish,  from  the  data,  signal  models  and  to  meas¬ 
ure  their  parameters,  as  contrasted  to  the  "direct  modeling"  waveform  methods  that  attempt  to 
fit  the  data  to  some  a  priori  model  without  any  attempt  to  verify  its  validity. 

The  purpose  of  the  deconvolution  methods  is  to  transform  the  complex  wavelets  consti¬ 
tuted  from  the  attenuation  operator,  instrument,  and,  in  some  methods,  the  estimate  of  the 
source  wavelet  and/or  the  site  function,  into  something  resembling  a  delta  function  which  can 
directly  be  recognized  in  the  output  by  visual  inspection  (Figure  1).  Preferably,  we  also  strive 
to  decrease  the  length  of  the  wavelet,  i.e.  increase  the  time  resolution.  However,  in  practice  it 
is  impossible  to  obtain  exact  delta  functions  by  deconvolution  because  of  the  bandwidth  limita¬ 
tions  in  the  signals,  the  fact  that  large  explosions  often  do  not  contain  enough  high  frequency 
energy,  and  the  ever  present  background  noise,  and  we  must  be  satisfied  with  approximations 
to  the  delta  function  such  as  broader  pulses  resembling  sin(x)/x  functions.  The  resulting  pulse 
is  often  broader  that  one  would  desire  and  small  pP  delay  times  may  not  always  be  resolvable. 
However,  such  limitations  are  fundamental  in  nature,  inherent  in  the  uncertainty  principle,  and 
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Figure  1.  Band-limited  deconvolution  results  for  synthetic  data.  Input  spike  sequences  for  a 
range  of  pP  times  are  shown  on  the  left  and  the  deconvolution  output  is  shown  on  the  right. 
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constitute  a  limit  to  the  applicability  of  the  seismic  method  which  cannot  be  overcome.  If  the 
pP  delay  time  is  smaller  that  die  width  of  the  resolution  kernel,  deconvolution  will  yield  pro¬ 
gressively  smaller  amplitudes  in  the  deconvolved  P+pP  waveforms  with  decreasing  pP  delay 
times  without  seeing  any  decrease  between  the  positive  and  negative  peak  distances  commonly 
used  to  estimate  pP  delay  times.  This  problem  is  illustrated  in  Figure  2.  Of  course,  for  actual 
nuclear  explosions  there  is  a  lower  limit  to  the  pP  delay  time,  assuming  everything  to  be  linear. 
This  is  imposed  by  the  minimum  depth  of  containment  in  the  rocks  in  question.  Moreover,  at 
very  shallow  depth  the  pP  phase  is  expected  to  disappear,  such  for  cratering  explosions.  If  the 
pP  arrival  can  be  resolved  in  the  deconvolutions  the  most  logical  approach  is  to  use  some 
measure  of  the  energy  in  the  first  positive  pulse  of  the  deconvolved  record  to  estimate  yield 
without  ever  estimating  mb.  This  approach  was  followed  by  Douglas  (1986)  and  should  be 
preferable  to  estimating  yields  from  mb.  The  advantage  of  the  maximum  likelihood  multichan¬ 
nel  deconvolution  method  discussed  in  this  paper  is  that  it  assumes  little  about  the  nature  of  the 
source  and  site  factors,  and  it  exploits  the  factorable  nature  of  P  wave  spectra  which  is  verified 
subsequently  in  each  case  by  the  reconstruction  of  the  original  data. 
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Figure  2.  (Top)  The  spectral  band  of  the  P  and  pP  wavelet  used  in  the  simulations.  (Bottom) 
Estimated  pP  delay  times  and  peak-to-peak  amplitudes  for  band  limited  deconvolutions  as  func¬ 
tions  of  true  pP  time  delay  (horizontal  axis). 
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DESCRIPTION  OF  THE  DECONVOLUTION  METHOD 

A  major  advance  in  the  understanding  of  the  structure  of  teleseismic  P  waves  was  the 
finding  by  Filson  and  Frasier  (1972)  that  P  wave  amplitude  spectra  from  a  set  of  events 
observed  at  a  narrow  range  of  azimuths  at  multiple  array  sites  may  be  factored  as 

Y ij(co)  =  Rj(to)  Sj(oo)  ,  (2) 

where  Y^to)  are  the  P  wave  spectra,  Sj(co)  are  source  factors,  and  Rj((0)  are  factors  appropriate 
to  the  receiver  sites.  The  validity  of  this  equation  has  been  verified  for  a  number  of  seismic 
arrays  as  well  as  more  widely  separated  receivers  (Lundquist  et  al  1980)  and  the  idea  was  also 
exploited  for  deriving  "relative  receiver  functions",  which  are  time  domain  representations  of 
intersite  transfer  functions  (Lundquist  et  al  1980)  including  Q  effects.  This  model  simply 
states  the  factorability  of  the  spectral  matrix  without  assuming  any  model  such  as  P+pP  as 
mentioned  above.  The  source  factors  should  include  all  features  that  associate  with  a  given 
source,  such  as  source  pulse  shape  (if  not  factored  out),  the  P,  pP,  and  spall  arrivals,  and  scat¬ 
tered  P  waves  generated  near  the  source,  such  as  Rayleigh-to-P  conversions  (Greenfield  1971). 
The  site  factors  include  peculiarities  due  to  the  geology  near  the  particular  sensor  in  the  array, 
such  as  P-to-Rayleigh  conversions  on  local  topography,  i.e.  anything  that  determines  the 
features  in  P  waveforms  peculiar  to  a  sensor  location.  The  two  factors  together,  as  we  shall 
show  below,  define  most  of  the  details  in  the  observed  P  waveforms  that  can  be  associated 
with  either  the  given  source  or  recording  site  provided  that  the  factorability  condition  as 
defined  by  equation  (2)  is  satisfied  for  the  given  data  set. 

In  this  paper  we  talk  of  "receiver  factors"  and  "source  factors"  in  order  to  make  the  dis¬ 
tinction  from  the  related  concept,  the  "relative  receiver  functions",  which  also  incorporate  some 
additional  assumptions  (Hart  et  al  1979),  such  as  the  maximum  "spikiness"  condition  (Wiggins 
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1978),  which  is  not  required  by  the  physics  of  the  problem.  We  separate  the  Q  operator,  the 
instrument,  and  the  source  wavelet  from  the  "source  factor"  and  combine  these  in  an  assumed 
known  source  wavelet,  Aj.  Thus,  Sj(co)  =  Aj(co)Xj(co)  where  A/co)  is  the  spectrum  of  the  com¬ 
bined  source  wavelet,  instrument,  and  Q.  This  way  the  "source  factor",  X/co),  is  a  Fourier 
transform  of  a  spike  sequence  associated  with  the  radiation  from  a  source  and  the  "receiver 
factor",  R;((0),  is  the  Fourier  transform  of  the  site  transfer  function.  Since  the  receiver  factors 
strongly  depend  on  the  azimuths  and  distances  of  the  sources,  they  may  be  assumed  to  be  con¬ 
stant  for  a  given  sensor  only  for  limited  source  regions. 

If  we  have  N  receivers  (sites)  and  M  sources,  we  may  describe  (  N  x  M  )  P  wave  spectra 
(or  time  series)  in  terms  of  N  +  M  factors;  there  is  thus  a  considerable  degree  of  redundancy 
in  the  problems  with  N  and  M  sufficiently  large.  This  helps  us  to  separate  the  R/co)  from  the 
Xj(co).  This  cannot  be  achieved  uniquely,  however,  unless  something  is  known  a  priori  about 
the  properties  of  each.  Otherwise,  common  spectral  factors  can  be  assigned  to  either  the 
source  or  receiver  factors,  or  canceling,  reciprocal  factors  may  be  assigned  to  both  the  receiver 
and  source  factors  as  indicated  below: 

Yij(G))  =[Ri(o)  F(co)]  [Fl(0))  X/co)]  A/co)  (3) 

This  is  the  most  dangerous  aspect  in  any  simultaneous  determination  of  source  and  site 
factors.  Fortunately,  if  we  estimate  the  source  factors  first  by  averaging  over  receiver  sites  we 
may  expect  that  in  the  first  step  the  site  effects  will  be  approximately  averaged  out.  This  is 
based  on  the,  quite  reasonable,  assumption  that  the  site  factors  will  tend  to  be  quite  random 
and  dissimilar  at  the  various  array  sensors.  In  the  following  steps  the  factors  are  refined 
without  too  much  change  and  thus  the  chances  of  introducing  major  distortion  into  the 
estimated  source  and  site  spike  sequences  is  thus  minimized. 
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Since  all  P  wave  signals  are  embedded  in  the  noise  background  it  is  desirable  that  the 
deconvolution-factoring  process  is  optimized  in  some  sense  with  respect  to  the  noise  i.e.,  we 
are  trying  to  find  some  compromise  between  the  resolution  and  the  S/N  ratio.  Thus,  although 
the  original  derivation  of  Shumway  (1984)  includes  weighting  with  respect  to  noise  spectra,  in 
this  work  we  do  not  maximize  the  S/N  ratio  in  order  to  increase  the  time  resolution  of  the 
deconvolution. 

The  detailed  derivation  of  the  maximum  likelihood  multichannel  deconvolution  method  is 
given  in  Der  et  al  (1983)  and  Shumway  and  Der  (1985).  Only  the  key  formulas  are  repeated 
here: 


where 


X  Aj  Ri  Y„ 


S I  Aj  |  2 1  Rs  |  2  +  Oj 


where  0J  =  — -  , 
P*j 


fc  I  Aj  |  2 1  Rs  |  2  +  0j 


and 


XAj  Xj  Yy 

_ i _ 

X I  Aj  |  2  [ I  Xj  |  2  +  d2 ] 
i 


Yjj  data 

Aj  source  wavelet  spectrum 
Rt  receiver  spike  response 
Xj  source  spike  sequence 
0  adjustable  N/S  ratio 


(4) 


(5) 


(6) 
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P„  noise  power  spectrum 

Pz  initial  signal  power 

0j2  effective  noise  power  in  second  step 

and  an  overbar  denotes  the  complex  conjugate.  The  i’s  correspond  to  the  individual  receiver 
sites  and  the  j’s  to  the  individual  sources.  Equations  4  through  6  are  very  similar  to  those 
presented  by  Oldenburg  (1981)  and  Oldenburg  et  al  (1981),  and  are  essentially  multichannel 
versions  of  the  same.  The  term  6  corresponding  to  the  signal-to-noise  ratio  performs  the  same 
function  as  in  Oldenburg’s  (1981)  formulas.  It  is  apparent  from  Equation  (4)  that  the  pro¬ 
cedure  is  also  similar,  in  computational  approach,  to  the  iterative  beaming  method  (Blandford 
et  al  1973)  in  which  the  previous  estimates  of  the  two  interfering  signals  are  iteratively 
estimated  in  each  step. 

The  deconvolution  is  done  in  the  frequency  domain.  The  first  step  is  to  factor  out  the 
"known"  spectra,  A,((o),  consisting  of  t*.  the  instrument  response,  and  the  source  function  (a 
von  Seggem  and  Blandford  (1972)  explosion  source  wavelet).  Since  we  are  working  with 
small  source  and  receiver  arrays,  we  assume  that  the  attenuation  parameter  t*P  is  the  same  for 
all  source-receiver  pairs  and  that  all  instruments  have  the  same  response.  At  the  start  of  the 
iteration  process  the  source  and  site  functions  are  each  initialized  to  an  impulse  in  the  time 
domain. 

The  iteration  starts  with  the  estimation  of  the  source  factors  by  deconvolving  the  source 
pulse  and  stacking  these  deconvolutions  over  the  available  sites.  This  is  intuitively  the  most 
reasonable  first  step  because  we  expect  the  site-generated  pulse  sequences  to  be  more  random 
than  the  source-related  ones  and  thus  they  may  be  easier  to  average  out  by  stacking.  The  flow 
diagram  of  the  iterations  is  depicted  in  Figure  3.  After  several  iterations,  each  source  factor 
will  contain  the  contribution  appropriate  to  the  particular  source  at  all  receivers  including  the 
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Figure  3.  Flow  diagram  of  the  multichannel  deconvolution  process.  In  the  center,  the  second 
and  third  boxes  from  the  top  correspond  to  Equations  (4)  and  (6)  in  the  text,  respectively. 
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source  pulse,  pP,  and  spall  (in  the  case  of  explosion  sources)  arrivals,  and  any  scattering  con¬ 
tributions  that  may  be  associated  with  a  particular  source. 

Each  site  factor  will  contain  the  contribution  of  that  site  to  the  seismogram,  mostly  the 
random  scattering  effects  attributable  to  the  structure  near  the  receiver.  The  site  factors  are  of 
high  frequency  in  character  to  compensate  for  the  partial  whitening  of  the  spectra  in  estimating 
the  source  time  functions  (Figures  4  and  5).  In  order  to  see  what  effects  these  have  on  the 
coda  one  has  to  extract  the  low  frequency  part  of  them  by  band-pass  filtering  (Figure  6).  The 
approximately  Sinc(x)  shaped  resolution  kernels  characterizing  the  attainable  resolution  with 
the  given  bandwidth  are  also  byproducts  of  the  iterative  process  (see  Shumway  and  Der  1985). 

As  a  final  test,  the  results  are  checked  by  reconstructing  each  input  trace  from  the  convo¬ 
lution  of  the  appropriate  source  and  site  terms  and  the  wavelet  consisting  of  the  von  Seggem 
and  Blandford  source  pulse,  t*,  and  instrument  response.  In  most  cases  the  reconstructions  and 
the  original  traces  are  very  similar  (Figures  7  and  8).  The  differences  between  the  original  and 
reconstructed  seismograms  can  be  inferred  as  a  measure  of  the  scattering  unique  to  each  indivi¬ 
dual  path  i.e.,  energy  unaccounted  for  by  either  site  or  source  factors,  and  thus  missing  in  the 
reconstructions. 

This  multichannel  deconvolution  approach  for  estimating  source  and  site  region  charac¬ 
teristics  has  several  advantages.  First,  the  method  utilizes  both  amplitude  and  phase  informa¬ 
tion  in  the  frequency  domain  calculations.  Furthermore,  no  a  priori  assumptions  need  to  be 
made  about  the  complexity  of  either  the  source  or  the  site  spike  series,  e.g.  the  presence  of  pP 
does  not  need  to  be  presupposed  by  providing  an  initial  guess  of  the  pP  delay  time  and  ampli¬ 
tude. 

These  iterative  formulas  yield  progressively  improved  estimates  of  the  source  and  receiver 
factors  and  their  equivalent  time  functions  as  evidenced  by  the  gradually  improving  quality  of 
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Figure  4.  Estimated  site  time  functions  (wrapped  around  time  zero)  at  three  sensors  of  YKA 
for  two  test  sites. 
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Figure  5.  Estimated  site  time  functions  at  three  sensors  of  EKA  for  three  test  sites. 
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Figure  6.  Original  (top)  and  band  pass  filtered  (bottom)  site  impulse  responses  at  EKA  for 
two  sensors  obtained  by  using  a  von  Seggem-Blandford  wavelet  as  a  convolution  operator  to 
extract  the  low-to-intermediate  frequency  contribution  to  the  site  response  function.  The  top 
sensor,  Rl,  has  a  considerable  amount  approximately  1  Hz  coda  due  to  the  site  response. 
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Figure  7a.  Examples  of  trace  reconstructions  (bottom)  and  original  traces  (top)  for  the  NTS 
event  Inlet  at  NORSAR  for  three  sites.  On  the  left  side  of  each  pair  of  traces  is  the  name  of 
the  sensor  and  the  correlation  coefficient  between  the  original  and  reconstructed  traces. 
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Figure  7b.  Examples  of  trace  reconstructions  (bottom)  and  original  traces  (top)  for  the  NTS 
event  Stilton  at  NORSAR  for  three  sites.  On  the  left  side  of  each  pair  of  traces  is  the  name  of 
the  sensor  and  the  correlation  coefficient  between  the  original  and  reconstructed  traces. 
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Figure  8a.  Examples  of  trace  reconstructions  (bottom)  and  original  traces  (top)  for  the  27 
April,  1975  E-EKTS  event  at  EKA  for  three  sites.  On  the  left  side  of  each  pair  of  traces  is  the 
name  of  the  sensor  and  the  correlation  coefficient  between  the  original  and  reconstructed  traces. 
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Figure  8b.  Examples  of  trace  reconstructions  (bottom)  and  original  traces  (top)  for  the  15 
January,  1965  cratering  E-EKTS  event  at  EKA  for  three  sites.  On  the  left  side  of  each  pair  of 
traces  is  the  name  of  the  sensor  and  the  correlation  coefficient  between  the  original  and  recon¬ 
structed  traces. 
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trace  reconstructions.  We  have  done  numerous  theoretical  simulations  with  known  inputs  and 
the  method  converged  in  all  cases  where  the  site  responses  were  set  to  be  notably  different 
from  the  initial  impulse  at  the  beginning  of  the  site  related  pulse  sequence.  We  have  found 
that  four  iterations  were  sufficient  to  stabilize  the  shapes  of  the  source  and  site  related  pulse 
sequences. 

Another  approach  is  to  include  only  the  t*  and  the  instrument  response,  and  not  the  von 
Seggem  and  Blandford  source  wavelet,  into  the  "known”  factor  to  be  divided  out  in  each  step 
of  the  iteration  (Douglas  1986).  This  way  the  source  pulses  may  directly  be  seen  in  the  source 
deconvolutions,  together  with  the  pulse  width  and  the  "overshoot".  This  approach  thus  has 
some  advantages  over  the  one  that  includes  an  estimate  of  the  explosion  source  wavelet, 
because  the  source  pulse  shapes  may  be  directly  seen,  but  it  is  harder  to  estimate  small  closely 
spaced  P  and  pP  delay  times  from  such  results.  This  approach  is  also  not  subject  to  errors  due 
to  incorrect  a  priori  assumptions  about  the  shapes  of  the  explosion  time  functions,  as  the 
explosion  source  parameters  (k,  B)  are  not  very  well  known  in  advance  (von  Seggem  and 
Blandford  1972,  Mueller  and  Murphy  1970,  Lay  et  al  1984b);  on  the  other  hand,  it  appears  the 
we  know  the  shapes  of  the  time  domain  attenuation  and  instrument  response  operators  ade¬ 
quately  enough  to  obtain  reliable  deconvolutions.  An  additional  argument  for  the  omission  of 
the  explosion  source  pulse  from  the  deconvolution  process  is  that  inclusion  of  such  a  pulse 
makes  sense  only  for  finding  the  pP  reflection  (that  has  supposedly  the  same  shape  as  the 
direct  P);  other  arrivals  (spall  and  arrivals  in  the  coda)  cannot  be  expected  to  have  the  same 
shape  as  the  direct  P  pulse.  If  the  deconvolution  is  done  without  removing  the  explosion  time 
function,  we  essentially  get  a  time  history  of  teleseismic  displacement  smoothed  by  the  approx¬ 
imately  sin(x)/x  shaped  resolution  kernel.  If  the  results  do  not  resemble  the  simple  waveforms 
expected  by  a  P+pP  model  (and  they  often  do  not)  then  the  conclusion  is  inescapable  that  such 
models  are  not  valid.  In  our  work  we  often  performed  the  deconvolutions  in  both  ways,  and 
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Teledyne  Geotech 


23 


October  1986 


ANNUAL  REPORT 


Multichannel  Deconvolution  of  P  Waves 


TGAL-86-6 


GENERAL  DATA  PROCESSING  PROCEDURES 

AH  traces  were  lined  up  in  time  on  the  first  arrival.  Since  we  are  not  interested  in  the 
absolute  amplitudes  in  the  context  of  estimating  relative  pP  times  and  amplitudes,  and  because 
we  want  to  give  roughly  equal  weight  to  all  the  traces  to  ensure  effective  factoring,  we  multi¬ 
plied  each  trace  with  two  factors,  one  appropriate  to  the  event,  the  other  one  to  the  site,  which 
were  designed  to  approximately  equalize  all  peak  amplitudes.  This  approach  is  consistent  with 
the  multiplicative  structure  of  the  spectra,  and  the  exact  values  of  the  equalizing  factors  do  not 
matter.  After  processing,  the  initial  normalization  factors  can  be  divided  out  to  estimate  the 
proper  relative  amplitudes  of  the  deconvolved  traces.  The  equalizing  of  traces  conflicts  some¬ 
what  with  the  optimization  of  the  S/N  ratios  in  the  output,  but  in  deconvolving  the  data 
improving  the  S/N  ratio  is  subordinated  to  the  goal  of  improving  the  time  resolution. 

As  mentioned  above,  the  iteration  is  started  with  the  initialization  of  the  time  domain 
representations  of  the  source  and  site  factors  to  delta  functions  and  the  source  factors  are 
estimated  first  by  summing  the  deconvolved  results  over  the  sites.  During  the  subsequent  itera¬ 
tions  both  the  source  and  site  factors  will  be  progressively  refined  by  factoring  the  spectral 
matrix.  The  first  estimates  of  source  functions  are  similar  to  the  stacked  deconvolutions 
obtained  by  Marshall  (1972)  and  Douglas  et  al  (1986).  Instead  of  the  frequency  domain 
approach  used  here  they  used  Wiener  filters  designed  in  the  time  domain  for  removing  t* 
effects  and  the  instrument  response.  Despite  these  differences  in  the  procedures  their  results 
usually  resemble  ours  closely. 

In  order  to  ensure  the  desired  deconvolution  effect  we  do  not  divide  by  the  noise  spectra, 
which  would  have  resulted,  especially  for  events  which  have  the  maxima  in  S/N  ratio  below  1 
Hz,  in  no  visible  deconvolution  at  all.  Instead,  we  have  assigned  a  unit  weight  to  all  frequen 
cies  which  had  S/N  ratios  above  I  and  excluded  the  frequency  ranges  with  poor  S/N  ratios 
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from  the  iterations.  Again,  this  degrades  the  S/N  in  the  output,  but  as  we  stated  above,  our 
main  goal  is  increased  resolution. 

The  results  require  some  interpretation:  one  must  keep  in  mind  that  the  deconvolution  of 
the  "known”  source  pulses  (consisting  of  t*.  instrument,  and  explosion  scaling  effects)  will 
produce  a  one  sided  pulse,  ideally  symmetric  in  time  and  having  some  side  lobes  because  of 
the  limitations  in  the  bandwidth  of  the  deconvolution  and  the  noncausal  nature  of  the  process. 
In  any  case,  the  results  are  easier  to  interpret  than  the  original  seismograms  which  are  complex 
superpositions  of  complex  wavelets.  The  resolution  is  typically  more  than  0.15  sec  which  is  the 
width  at  half  maximum  of  the  resolution  kernels.  The  deconvolved  site  terms  are  byproducts 
of  the  deconvolution  process  in  the  sense  that  we  are  mostly  interested  in  the  source  time  func¬ 
tions  in  this  study.  As  we  pointed  out  above,  since  the  optimized  deconvolutions  for  source 
terms  cannot  have  white  spectra  because  of  noise  constraints,  the  "site”  factors  contain  a  lot  of 
high  frequency,  much  of  which  is  just  noise,  and  only  the  low  frequency  pan  of  which  is  phy¬ 
sically  meaningful  The  site  time  functions  also  appear  "wrapped  around"  around  time=0  sec 
because  of  the  way  the  factoring  was  set  up.  Some  precursors  are  visible  because  the  process 
does  not  have  a  causality  constraint,  and  partly  because  of  the  wraparound  effect.  The  subse¬ 
quent  reconstructions  of  the  original  traces  demonstrate,  however,  that  noncausality  is  not  a 
major  problem  since  most  reconstructed  traces  do  not  have  major  precursors. 

In  all  our  deconvolution  work  we  have  used  only  complete  matrices  of  event  by  site  of  P 
vjve  data  For  a  given  array  often  several  such  data  sets  can  be  found  that  have  only  a  lim¬ 
ited  subset  of  common  events  and  sensors  sites.  It  is  reassuring,  however,  that  for  the  com¬ 
mon  events  the  source  time  function  estimates  are  commonly  quite  similar  regardless  of  the 
particular  data  set  This  indicates  the  robustness  of  the  method. 
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The  deconvolutions  presented  in  the  following  sections  have  all  utilized  the  best  constant 
t*  operator  for  a  given  path  as  opposed  to  a  frequency  dependent  t*  operator,  t*(f).  This  has 
been  shown  previously  (i.e.  Der  and  Lees  1985)  to  be  a  reasonable  way  to  treat  narrow  band 
data  when  t*(f)  is  slowly  varying.  For  t*  and  t*(f)  models  that  are  appropriate  for  the  given 
path,  the  deconvolved  source  time  functions  show  only  minor  differences  in  waveform  shape. 
The  absolute  amplitudes  of  the  source  functions  do  change,  and  this  may  be  important  in  yield 
estimation  (to  be  discussed  later). 

The  quality  of  the  reconstructions  (and  the  validity  of  the  spectral  factorization  assump¬ 
tion)  is  evaluated  in  each  case  by  computing  the  correlation  coefficients  between  the  recon¬ 
structed  and  original  waveforms.  In  the  rest  of  this  report  we  shall  only  show  a  small  fraction 
of  such  results  and  only  when  it  can  be  justified  in  making  an  important  point.  The  work  has 
produced  massive  amounts  of  output  that  cannot  be  shown  in  a  report  of  this  size.  All  results 
are  available  at  request. 
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RESULTS  OF  THE  DATA  ANALYSES 

Thus  far  in  this  project  we  have  processed  the  following  data  sets:  C-EKTS  explosions 
plus  one  earthquake  as  recorded  at  the  AWRE  arrays  EKA,  GBA,  WRA,  and  YKA,  E-EKTS 
explosions  recorded  at  EKA,  GBA,  WRA,  and  YKA,  Novaya  Zemlya  events  recorded  at  EKA, 
GBA,  and  YKA,  Yucca  Flats  events  plus  Piledriver  recorded  at  EKA,  GBA,  and  YKA,  Pahute 
Mesa  events  including  Faultless  recorded  at  EKA,  GBA,  and  WRA,  and  French  Sahara  explo¬ 
sions  recorded  at  EKA  and  YKA.  We  also  present  preliminary  results  for  E-EKTS  explosions 
deconvolved  at  the  3-component  RSTN  stations.  The  events  in  each  data  set  are  described  in 
Tables  la  to  If. 

For  describing  the  test  sites  within  the  Soviet  East  Kazakhstan  test  site,  we  are  using  the 
descriptors  E-EKTS,  C-EKTS,  and  NW-EKTS  as  defined  in  Figure  9  for  the  eastern,  central, 
and  northwestern  test  sites  within  East  Kazakhstan.  These  test  sites  are  described  and  named 
on  pages  182  and  183  of  Dahlman  and  Israelson  (1977)  and  in  Rodean  (1979).  Dahlman  and 
Israelson  describe  the  E-EKTS,  C-EKTS,  and  NW-EKTS  test  sites  as  a  marshy  lowlands,  a 
mountainous  region,  and  a  plain,  respectively. 

For  nearly  all  of  the  deconvolutions  we  have  used  data  recorded  at  the  AWRE  arrays 
EKA,  GBA,  WRA,  and  YKA.  These  arrays  are  described  by  Mowat  and  Burch  (1977).  Table 
2  lists  some  of  the  parameters  of  the  arrays  and  layouts  of  the  arrays  are  shown  in  Figures  10 
a-d. 


Results  for  C-EKTS  Explosions 

A  map  of  the  locations  of  the  C-EKTS  explosions  is  shown  in  Figure  11.  The  data  was 
deconvolved  at  each  array  using  the  options  of  removing  and  not  removing  the  estimated  von 
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Table  la 


Central  EKTS  Events  Used  in  Deconvolutions 


Event 

Array 

Location 

EKA 

YKA  WRA 

GBA 

Origin  Time 

Lat(°N) 

Lon(°E) 

mb 

670226 

X 

03:57:57.4 

49.750:78.125 

6.0 

711230 

X 

06:20:57.8 

49.772 

78.093 

5.8 

770329 

X 

X 

X 

03:56:57.6 

49.790 

78.086 

5.4 

770730 

X 

X 

X 

X 

01:56:57.7 

49.759 

78.097 

5.1 

780326 

X 

X 

X 

03:56:57.6 

49.768 

78.044 

5.6 

780422 

X 

X 

X 

X 

03:06:57.6 

49.761 

78.186 

5.3 

781031 

X 

X 

X 

X 

04:16:57.8 

49.806 

78.143 

5.2 

790506 

X 

X 

X 

X 

03:16:57.7 

49.774 

78.049 

5.2 

790531 

X 

X 

X 

X 

05:54:57.7 

49.835 

78.127 

5.2 

791018 

X 

X 

04:16:57.7 

49.837 

78.148 

5.2 

800522 

X 

X 

X 

03:56:57.7 

49.784 

78.082 

5.5 

800731 

X 

X 

X 

X 

03:32:57.7 

49.807 

78.148 

5.3 

Location  information  from  Marshall  et  al  (1984). 


Table  lb 


Eastern  EKTS  Events  Used  in  Deconvolutions 


Event 

Array 

Location 

EKA 

YKA  WRA 

GBA 

Origin  Time 

Lat(°N) 

Lon(°E) 

"Y 

650115 

X 

05:59:58.4 

49.940 

79.010 

5.8 

720210 

X 

05:02:57.5 

50.014 

78.878 

5.4 

721210 

X 

04:27:07.3 

50.001 

78.973 

6.0 

731214 

X 

07:46:57.2 

50.044 

78.987 

5.8 

741016 

X 

X 

X 

06:32:57.6 

49.979 

78.898 

5.5 

750630 

X 

03:26:57.6 

50.004 

78.957 

4.6 

751029 

X 

X 

04:46:57.3 

49.946 

78.878 

5.7 

751225 

X 

X 

X 

05:16:57.2 

50.044 

78.814 

5.7 

760421 

X 

X 

X 

05:02:57.2 

49.890 

78.827 

5.3 

760609 

X 

X 

03:02:57.2 

49.989 

79.022 

5.1 

760704 

X 

02:56:57.5 

49.909 

78.911 

5.8 

761123 

X 

05:02:57.3 

50.008 

78.963 

5.9 

761207 

X 

04:56:57.4 

49.922 

78.846 

5.9 

770629 

X 

X 

X 

03:06:57.8 

50.006 

78.869 

5.2 

770905 

X 

X 

X 

03:02:57.3 

50.035 

78.921 

5.7 

771130 

X 

04:06:57.4 

49.958 

78.885 

6.0 

780611 

X 

X 

X 

02:56:57.6 

49.898 

78.797 

5.9 

780705 

X 

X 

02:46:57.5 

49.887 

78.871 

5.8 

780829 

X 

X 

X 

02:37:06.3 

50.000 

78.978 

5.9 

780915 

X 

X 

02:36:57.4 

49.916 

78.879 

6.0 

781104 

X 

X 

X 

05:05:57.3 

50.034 

78.943 

5.6 

781129 

X 

04:33:02.5 

49.949 

78.798 

6.0 

790623 

X 

02:56:57.5 

49.903 

78.855 

6.2 

790707 

X 

03:46:57.3 

50.026 

78.991 

5.8 

790804 

X 

03:56:57.1 

49.894 

78.904 

6.2 

791028 

X 

03:16:56.9 

49.973 

78.997 

6.0 

791202 

X 

04:36:57.5 

49.891 

78.796 

6.0 

791223 

X 

04:56:57.4 

49.916 

78.755 

6.2 

800629 

X 

02:32:57.7 

49.939 

78.815 

5.7 

801227 

X 

04:09:08.1 

50.057 

78.981 

5.9 

810422 

X 

01:17:11.3 

49.885 

78.810 

6.1 

811018 

X 

03:57:02.6 

49.923 

78.859 

6.1 

811129 

X 

03:35:08.6 

49.887 

78.860 

5.7 

820831 

X 

01:31:00.7 

49.924 

78.761 

5.4 

821226 

X 

03:35:14.2 

50.071 

78.988 

5.7 

Location  information  from  Marshall  et  al  (1984). 
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Novaya  Zemlya  Events  Used  in  Deconvolutions 


Event 

Array 

Location 

EKA 

YKA  WRA 

GBA 

Origin  Time 

Lat(°N) 

Lon(°E) 

mb 

730927 

X 

X 

06:59:58.5 

70.756 

53.754 

5.84 

760927 

X 

02:59:57.7 

73.360 

54.880 

5.77 

761020 

X 

X 

07:59:58.1 

73.399 

54.835 

4.89 

770901 

X 

X 

02:59:58.0 

73.339 

54.626 

5.71 

780810 

X 

07:59:58.0 

73.293 

54.885 

6.04 

780927 

X 

X 

02:04:58.6 

73.350 

54.677 

5.68 

790924 

X 

X 

03:29:58.8 

73.346 

54.679 

5.80 

791018 

X 

X 

07:09:58.8 

73.318 

54.821 

5.85 

801011 

X 

X 

X 

07:09:57.5 

73.335 

54.938 

5.80 

811001 

X 

X 

12:14:57.3 

73.308 

54.817 

5.91 

821011 

X 

X 

X 

07:14:58.7 

73.348 

54.601 

5.52 

830818 

X 

X 

16:09:58.9 

73.358 

55.974 

5.84 

830925 

X 

X 

13:09:58.2 

73.326 

54.564 

5.71 

Location  information  from  Marshall  et  al  (1986). 
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Table  Id 


French  Sahara  Events  Used  in  Deconvolutions 

Event  Array 

EKA  YKA  WRA  GBA 
Emeraude  X  X 

Grenat  X  X 

Rubis  X 

Saphir  X  X 


Table  le 


Pahute  Mesa  Events  Used  in  Deconvolutions 


Event 

Array 

Location 

EKA 

YKA  WRA 

GBA 

Date 

Origin  Time 

Lat(°N) 

Lon(°W) 

BOXCAR 

X 

680426 

15:00:00.0 

37.295 

116.456 

CAMEMBERT 

X 

X 

750626 

12:30:00.2 

37.279 

116.369 

CHESHIRE 

X 

760214 

11:30:00.2 

37.243 

116.420 

ESTUARY 

X 

X 

760309 

14:00:00.1 

37.310 

116.364 

FAULTLESS 

X 

680119 

18:15:00.1 

38.63 

116.22 

FONTINA 

X 

760212 

14:45:00.2 

37.271 

116.488 

HANDLEY 

X 

700326 

19:00:00.2 

37.300 

116.534 

INLET 

X 

751120 

15:00:00.1 

37.225 

116.368 

KASSERI 

X 

X 

X 

751028 

14:30:00.2 

37.290 

116.412 

MAST 

X 

750619 

13:00:00.1 

37.350 

116.320 

MUENSTER 

X 

X 

X 

760103 

19:15:00.2 

37.297 

116.333 

PURSE 

X 

690507 

13:45:00.0 

37.283 

116.501 

TYBO 

X 

X 

750514 

14:00:00.2 

37.221 

116.474 

Locations  were  announced  except  for  FAULTLESS,  for  which  the  location  is  from  the  ISC. 


i 

\ 


32 


Table  If 


Yucca  Flats  Events  Used  in  Deconvolutioas 


Event  Array  Location 


EKA 

YKA 

WRA  GBA 

Date 

Origin  Time 

Lat(°N) 

Lon(°W) 

BULKHEAD 

X 

X 

770427 

15:00:00.1 

37.095 

116.028 

CABRILLO 

X 

750307 

15:00:00.2 

37.134 

116.084 

CHIBERTA 

X 

751220 

20:00:00.2 

37.128 

116.062 

CREWLINE 

X 

X 

770525 

17:00:00.1 

37.094 

116.045 

FARALLONES 

X 

X 

771214 

15:00:00.1 

37.136 

116.086 

HEARTS 

X 

790906 

15:00:00.1 

37.088 

116.053 

LOWBALL 

X 

X 

780712 

17:00:00.1 

37.079 

116.044 

QUARGEL 

X 

X 

781118 

19:00.00.2 

37.127 

116.084 

SCANTLING 

X 

770819 

17:55:00.1 

37.110 

116.055 

STRAIT 

X 

X 

760317 

14:45:00.1 

37.107 

116.052 

TOPGALLANT 

X 

X 

750228 

15:15:00.1 

37.106 

116.056 

PILEDRIVER 

X 

660602 

15:30:00.1 

37.227 

116.056 

Locations  were  announced. 
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a  known  fault.  There  are  three  main  testing  areas  within  the  Eastern  Kazakh  Tesat  Site 
(EKTS);  a  northwestern  diffuse  area,  a  central  area  covering  a  rugged,  mountainous  region  and 
an  eastern  area  clearly  separated  from  the  rest.  These  will  be  referred  to  by  the  abbreviations 
NW-EKTS,  C-EKTS  and  E-EKTS  in  the  rest  of  this  report. 


Table  2 


United  Kingdom  Sponsored  Arrays 

Maximum  Date  of 


Element 

Element 

Date 

Digital 

Station 

Code 

Location 

Spacing 

Spacing 

Operational 

Recording 

(km) 

(km) 

Eskdalemuir 

EKA 

55°19’59"N 

0.9 

9.8 

17  May  1962 

14  Nov  1983 

Scotland 

3°09‘33"W 

Yellowknife 

YKA 

62°29’34"N 

2.5 

22.5 

26  Nov  1962 

- 

Canada 

1 14°36’ 17" W 

Gauribidanur 

GBA 

13°36’ 15"N 

2.5 

32.0 

1  Feb  1966 

4  Mar  1979 

India 

77°26’ 10" E 

NVarramunga 

WRA 

19°56’39"S 

2.5 

26.3 

1  Mar  1966 

7  Jun  1977 

Australia  134°20’27"E 
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Figure  lOd.  Layout  of  the  AWRE  array  WRA. 
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Seggern  Blandford  (VSB)  explosion  source  wavelet.  Figures  12  through  15  show  the  source 
functions  at  the  various  arrays  with  the  estimated  VSB  removed  in  the  deconvolution  process. 
The  P  and  pP  phases  are  marked  on  each  trace.  It  is  apparent  that  only  a  few  of  the  events 
show  secondary  arrivals,  interpretable  as  pP,  with  negative  amplitudes  similar  in  size  to  the  ini¬ 
tial  P  wave  pulse.  Most  events  do  not  show  a  clear  pP.  It  may  be  argued  that  many  C-EKTS 
events  were  shallow  enough  for  preventing  the  development  of  pP  pulses,  possibly  by  non¬ 
linearity  of  deformation  near  the  surface.  An  alternative  explanation  may  be  the  complex 
topography  above  the  explosions. 

Another  interesting  feature  is  the  differences  in  the  deconvolved  waveforms  at  the 
different  arrays.  In  Figure  16  we  show  the  results  obtained  for  common  events  at  the  four 
AVVRE  arrays,  this  time  without  removing  the  VSB  source  pulse.  The  events  are  loosely 
arranged  by  the  locations  in  Figure  11.  For  events  recorded  at  a  given  array,  there  often 
appear  to  be  more  similarities  between  the  source  functions  of  events  that  were  spacially  close 
to  each  other,  though  we  have  not  studied  this  in  a  quantitative  manner.  The  event  of  March 
20,  1976  (760320)  is  an  earthquake  (see,  for  example,  Rodean  1979).  Variations  in  pP  and  sP 
amplitude  and  polarity  are  clearly  seen  at  the  different  azimuths  for  this  event.  Another 
interesting  event  (770329)  has  been  identified  as  a  double  event  (McLaughlin  et  al  1986c)  with 
time  shifts  appropriate  to  different  distances  relative  to  the  various  arrays.  For  this  event,  note 
that  the  second  event  has  a  clear  pP  arrival,  unlike  most  of  the  C-EKTS  events.  The  second 
event  is  an  explosion  which  occurred  at  the  NW-EKTS  test  site,  which  is  a  relatively  flat  and 
simple  geologic  region. 

With  these  events  aside  there  is  still  a  large  difference  among  the  deconvolved  waveforms 
at  the  different  arrays,  WRA  and  GBA  show  more  of  a  high  frequency  coda  than  the  other  two 
arrays,  (Some  of  the  complexity  at  WRA  may  be  due  to  its  distance  of  -85"  from  C-EKTS. 
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Figure  12.  Deconvolved  source  functions  for  C-EKTS  events  at  EKA.  The  estimated  von 
Seggem  and  Blandford  wavelet  has  been  removed  in  the  deconvolutions.  The  P  and  pP,  if 
present,  arrivals  are  marked  on  each  trace. 


Figure  13.  Deconvolved  source  functions  for  C-EKTS  events  at  GBA.  The  estimated  von 
Seggern  and  Blandford  wavelet  has  been  removed  in  the  deconvolutions.  The  P  and  pP,  if 
present,  arrivals  are  marked  on  each  trace. 
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Figure  14.  Deconvolved  source  functions  for  C-EKTS  events  at  WRA.  The  estimated  von 
Seggern  and  Blandford  wavelet  has  been  removed  in  the  deconvolutions.  The  P  and  pP,  if 
present,  arrivals  are  marked  on  each  trace. 
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Figure  15.  Deconvolved  source  functions  for  C-EKTS  events  at  YKA.  The  estimated  von 
Seggern  and  Blandford  wavelet  has  been  removed  in  the  deconvolutions.  The  I’  and  pi’,  if 
present,  arrivals  are  marked  on  each  trace. 
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Figure  16.  A  set  of  source  time  functions  for  the  C-EKTS  test  site  and  an  earthquake  (bot¬ 
tom)  obtained  at  AWRE  arrays.  VSB  wavelets  have  not  been  removed  in  these  deconvolu¬ 
tions.  Most  of  the  traces  do  not  show  a  clearly  identifiable  pP.  The  770329  event  consists  of 
two  nearly  coincident  events  at  two  test  sites  (C-EKTS  and  NW-EKTS)  and  shows  variable 
time  delays  at  the  various  arrays. 
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At  this  distance,  the  Jeffreys-Bullen  travel  time  curves  predict  that  PcP  will  arrive  ~3.8  sec 
after  the  P,  and  ir  his  case  perhaps  even  sooner  since  both  the  source  and  receiver  regions  are 
on  shields  which  tend  to  have  faster  than  average  velocities.)  There  are  two  working 
hypotheses  to  explain  this  phenomenon.  The  first  is  that  there  are  components  in  the  site 
responses  at  each  array  (multiplicative  factors  in  the  frequency  domain)  that  need  to  be 
removed  in  order  to  equalize  the  results.  The  other  is  that  the  time  functions  radiated  by  the 
various  events  are  not  the  same  in  all  directions,  i.e.  the  source  radiation  is  azimuthally  asym¬ 
metrical.  To  decide  which  of  these  hypotheses  is  true  we  have  selected  a  subset  of  events  and 
a  subset  of  sensors  at  each  array  and  jointly  deconvolved  the  data.  Joint  deconvolution  of 
multi-array  data  implies  that  the  assumptions  appropriate  to  the  first  hypothesis,  i.e.  that  the 
sources  are  azimuthally  symmetrical  and  that  the  differences  in  the  deconvolution  results  can 
be  removed  by  additional,  multiplicative  site  factors  applied  to  each  array,  are  incorrect.  Table 
3  lists  the  events  and  sensors  used  in  the  joint  deconvolution  and  Figure  17  shows  the  source 
functions  from  the  joint  deconvolutions  as  well  as  from  deconvolution  of  the  same  set  of 
events  separately  at  each  array.  The  source  functions  derived  from  the  joint  deconvolutions  are 
much  simpler  than  the  source  functions  from  the  deconvolutions  at  the  separate  arrays.  How¬ 
ever,  the  source  and  site  factors  obtained  from  the  joint  deconvolution  result  in  a  noticeable 
degradation  of  the  reconstructions  of  the  original  data  indicating  that  the  differences  cannot  be 
attributed  to  differences  in  the  site  functions.  This  is  shown  in  Figure  18  where  reconstruc¬ 
tions  using  the  source  and  site  factors  from  the  joint  deconvolution  give  correlation  coefficients 
of  only  0.4  to  0.8  when  compared  with  the  original  traces  while  reconstructions  using  source 
and  site  factors  from  deconvolution  of  the  same  events  at  a  single  array  give  correlation 
coefficients  of  greater  that  0.9.  This  leads  to  the  conclusion  that  for  the  C-EKTS  events,  the 
differences  of  the  deconvolution  results  at  the  various  arrays  must  be  attributed  to  azimuthal 
asymmetries  in  the  source  radiation  of  the  individual  events. 
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Table  3 


Events  and  Receivers  used  in  Joint  Deconvolution 
of  Central  EKTS  Events  at  AWRE  Arrays 


Sensors  Used  at  Each  Array 


Events 

EKA 

GBA 

WRA 

YKA 

770730 

R1 

R2 

R1 

R1 

781031 

R4 

R4 

R2 

R5 

790506 

R5 

B1 

R5 

R7 

780422 

R8 

B2 

R8 

R9 

790531 

B1 

B7 

B1 

B3 

800731 

B4 

B9 

B4 

B8 

48 


Figure  17.  Comparison  of  source  deconvolutions  (VSB  not  removed)  for  a  set  of  common 
events  at  each  of  the  AWRE  arrays  and  a  joint  deconvolution  using  six  sensors  at  each  array. 


Figure  18.  Comparison  of  some  original  traces  and  reconstructions  obtained  from  the  joint 
deconvolution  of  common  events  at  the  four  AWRE  arrays.  The  poor  quality  of  reconstruc¬ 
tions  (unlike  those  at  individual  arrays)  indicates  a  considerable  asymmetry  in  the  source  radia¬ 
tion  from  C-EKTS. 
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Results  for  E-EKTS  Explosions 

Various  subsets  of  the  E-EKTS  explosion  data  were  deconvolved  at  EKA,  GBA,  WRA, 
and  YKA.  Marshall  et  al  (1984)  classified  the  waveforms  from  E-EKTS  events  into  types  1 
(simpler)  and  2  (more  complex),  the  origin  locations  of  which  are  shown  on  the  map  in  Figure 
19.  Figure  20  is  a  map  of  the  locations  of  the  E-EKTS  events  deconvolved  in  this  study,  and 
Figures  21  to  24  show  the  source  functions  from  deconvolutions  with  the  VSB  source  pulse 
removed.  Most  events  show  strong  secondary  arrivals  of  negative  polarity  which  may  be  inter¬ 
preted  as  pP  following  the  initial  P  wave  arrival.  In  this  respect,  the  E-EKTS  results  are 
noticeably  different  from  those  obtained  for  the  C-EKTS  test  site.  One  exception  is  the  crater¬ 
ing  explosion  (event  650115  at  the  bottom  of  Figure  21)  which  does  not  have  a  clear  pP 
arrival.  This  event  was  also  shown  to  have  a  magnitude  bias  relative  to  other  (E-EKTS)  events 
used  for  comparison  (Der  et  al  1986a,b).  To  the  right  of  each  trace  in  Figures  21  to  24  is  a 
"1"  or  "2"  according  to  the  waveform  classification  of  Marshall  et  al  (1984)  for  of  E-EKTS 
events.  Class  "1"  events  show  simpler  P  +  pP  waveforms  than  do  Class  "2”  events,  and  the 
results  of  our  deconvolutions  are  generally  consistent  with  the  classifications  of  Marshall  et  al. 

One  may  ask  the  question  whether  the  differences  between  the  C-EKTS  and  E-EKTS 
source  functions  may  be  due  to  differences  between  the  site  responses  for  the  two  test  sites. 
Based  on  past  results  the  answer  must  be  negative.  The  E-EKTS  and  C-EKTS  test  sites  are 
not  distant  enough  to  cause  appreciable  differences  in  the  site  responses  at  common  arrays,  and 
it  is  even  less  likely  that  this  would  happen  for  three  arrays  at  the  same  time.  Moreover,  we 
have  been  successful  in  the  past  in  deconvolving  a  mix  of  C-EKTS  and  E-EKTS  events 
recorded  at  NORSAR  using  the  same  site  functions  for  both  test  sites,  and  the  trace  reconstruc¬ 
tions  were  of  excellent  quality  indicating  that  the  site  functions  were  practically  the  same  (Der 
et  al  1986a,b). 
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Figure  21  .  Source  time  function  estimates  obtained  for  E-EKTS  explosions  at  EKA.  A  VSB 
wavelet  has  been  removed  in  these  deconvolutions.  The  waveform  classification  assigned  by 
Marshall  et  al  (1984)  is  noted  to  the  right  of  each  trace.  Pulses  interpretabie  as  pP  are  typical 
of  E-EKTS  deconvolutions.  However,  note  that  the  pP  is  missing  for  the  cratering  event 
(6501 15). 
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Marshall,  Bache,  and  Liiwall  (1984) 
Waveform  Classification  of  Shagan  Events 
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Figure  21  Continued. 
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Marshall,  Bache,  and  Lilwall  (1984) 
Waveform  Classification  of  Shagan  Events 
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Marshall,  Bache,  and  LiKvall  (1984) 
Waveform  Classification  of  Shagan  Events 
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Figure  22.  Source  time  function  estimates  obtained  for  E-EKTS  explosions  at  GBA.  A  VSB 
wavelet  has  been  removed  in  these  deconvolutions.  The  waveform  classification  assigned  by 
Marshall  et  al  (1984)  is  noted  to  the  right  of  each  trace. 
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Figure  22  Continued. 
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Figure  23.  Source  time  function  estimates  obtained  for  E-EKTS  explosions  at  WRA.  A  VSB 
wavelet  has  been  removed  in  these  deconvolutions.  The  waveform  classification  assigned  by 
Marshal!  et  al  (1984)  is  noted  to  the  right  of  each  trace. 
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Figure  24.  Source  time  function  estimates  obtained  for  E-EKTS  explosions  at  YKA.  A  VSB 
wavelet  has  been  removed  in  these  deconvolutions.  The  waveform  classification  assigned  by 
Marshall  et  al  (1984)  is  noted  to  the  right  of  each  mace. 
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Results  for  Novaya  Zemlya  Events 

The  estimated  locations  of  the  13  Novaya  Zemlya  events,  analyzed  at  EKA,  GBA,  and 
YKA,  are  shown  in  Figure  25.  We  have  12  events  from  the  northern  Novaya  Zemlya  test  site 
and  one  event  from  the  southern  Novaya  Zemlya  test  site.  Figures  26  to  28  show  source  time 
functions  at  the  three  arrays  with  the  VSB  source  pulse  removed  and  Figure  29  is  a  composite 
of  deconvolutions  at  the  separate  arrays,  without  the  VSB  removed.  In  Figure  29,  the  events 
are  loosely  arranged  by  location  according  to  die  map  in  Figure  25.  Event  801011  has  been 
identified  by  Marshall  et  al  (1986)  as  a  double  event  and  in  our  deconvolutions  this  can  best  be 
seen  on  the  source  function  for  801011  as  deconvolved  at  EKA. 

The  most  striking  observation  about  the  waveforms  in  Figure  29  is  how  much  more 
complex  the  P  wave  and  P  wave  coda  are  at  YKA  than  at  EKA  and  GBA.  However,  even 
though  the  waveforms  at  YKA  are  very  complex,  they  do  reconstruct  very  well  as  shown  in 
Figure  30.  Greenfield  (1971)  has  previously  noted  the  complex  coda  for  Novaya  Zemlya 
events  recorded  at  YKA  and  McLaughlin  et  al  (1968c)  have  observed  that  stations  at  azimuths 
of  352  to  111°  from  Novaya  Zemlya  have  significantly  shorter  P  wave  coda  than  stations  at 
azimuths  of  150  to  349°.  While  the  deconvolved  sources  at  YKA  are  too  complex  to  identify 
a  pP  arrival,  a  pP  arrival  can  be  recognized  on  most  of  the  deconvolutions  of  recordings  at 
EKA  and  GBA.  Even  though  the  source  functions  at  EKA  and  GBA  are  much  simpler  than 
those  at  YKA,  they  consistently  show  late  arrivals  other  than  pP.  Some  of  the  later  arrivals  at 
EKA  may  be  due  to  triplications  in  the  P  wave  travel  time  curve  which  still  has  two  observ¬ 
able  branches  at  a  distance  of  29°  for  shield  models  such  as  that  of  King  and  Calcagnile 
(1976). 

It  is  also  interesting  to  examine  the  source  functions  for  the  one  event  from  southern 
Novaya  Zemlya  (730927),  This  event  has  simpler  source  functions  at  both  EKA  and  YKA 
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Figure  25.  Locations  of  Novaya  Zemlya  events  analyzed  in  this  report 
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Figure  26.  Source  time  function  estimates  obtained  for  Novaya  Zemlya  events  at  EKA.  A 
VSB  wavelet  has  been  removed  in  these  deconvolutions. 
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Figure  27.  Source  time  function  estimates  obtained  for  Novaya  Zemlya  events  at  GBA.  A 
VSB  wavelet  has  been  removed  in  these  deconvolutions. 


Figure  28.  Source  time  function  estimates  obtained  for  Novaya  Zemlya  events  at  YKA. 
VSB  wavelet  has  been  removed  in  these  deconvolutions. 
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Hj»un.*  29.  Source  lime  function  estimates  for  Novaya  Zemlya  events  at  three  AWRK  arrays 
A  VSB  wavelet  has  not  been  removed  in  these  deconvolutions.  The  source  time  functions  (as 
well  as  the  original  seismograms)  are  quite  complex  at  YKA  in  comparison  with  the  other 
arrays 
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Figure  30.  Examples  of  trace  reconstructions  compared  to  original  traces  at  selected  sites  of 
YKA.  The  quality  of  reconstruction  is  quite  good. 
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than  do  the  northern  Novaya  Zemlya  events.  One  possible  explanation  for  this  may  be  the 
difference  in  surface  topography  between  the  northern  and  southern  test  sites;  the  northern 
Novaya  Zemlya  test  site  is  very  mountainous  while  the  southern  test  site  is  in  a  much  Hatter 
area  (per  map  in  Marshall  et  al,  1986). 

Results  for  French  Sahara  Nuclear  Tests 

A  map  of  the  locations  of  the  main  events  at  this  test  site  is  shown  in  Figure  31.  Only  a 
few  nuclear  explosions  in  Sahara  were  available  with  acceptable  S/N  ratios  at  AWRE  arrays. 
The  Sahara  test  site  is  another  low  Q  area  (Der  et  al  1985),  thus  making  deconvolutions 
difficult.  The  deconvolutions  yielded  fairly  complex  time  functions  for  the  sources  which  vary' 
considerably  with  azimuth  (Figure  32)  but  with  good  reconstructions  in  all  directions.  We 
interpret  this  as  azimuthal  asymmetries  in  the  radiation  which,  in  the  view  of  the  complex 
topography  of  the  test  site  is  to  be  expected  (McLaughlin  et  al  1986b).  Due  to  the  limited 
bandwidth  of  the  records  we  did  not  attempt  to  factor  out  the  VSB  source  wavelet. 

Results  for  Pahute  Mesa  Events 

Twelve  Pahute  Mesa  events  plus  the  shot  Faultless  were  deconvolved  at  EKA,  GBA, 
WRA,  and  YKA.  A  map  of  the  relative  locations  of  the  events  is  shown  in  Figure  33  and 
deconvolved  source  functions,  with  the  VSB  source  pulse  removed,  are  shown  in  Figures  34  to 
36.  The  deconvolutions  of  the  Pahute  Mesa  events  show  a  great  deal  of  ringing  and  a  rela¬ 
tively  narrow  bandwidth  due  to  the  high  Q  in  the  NTS  region.  Reliable  pP  phases  cannot  be 
picked  for  most  of  the  events.  The  source  function  of  Faultless  deconvolved  at  EKA  (on  Fig¬ 
ure  34)  does  not  appear  significantly  different  from  the  deconvolved  source  functions  of  othet 
large  events  at  Pahute  Mesa. 
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Figure  31.  Locations  of  tests  at  the  French  Sahara  test  area. 
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re  32.  Deconvolution  results  for  source  time  estimates  for  all  the  available  data  from  the 
ch  Sahara  events  at  EKA  and  YKA.  VSB  wavelets  have  not  been  removed  in  these 
nvolutions. 
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Figure  34.  Source  time  function  estimates  for  thirteen  Pahute  Mesa  explosions  at  EKA.  VSB 
wavelets  have  been  removed  in  these  deconvolutions. 
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Figure  34  Continued. 
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Figure  35.  Source  time  function  estimates  for  five  Pahute  Mesa  explosions  at  GBA.  VSB 
wavelets  have  been  removed  in  these  deconvolutions. 
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Figure  36.  Source  time  function  estimates  for  two  Pahute  Mesa  explosions  at  WRA.  VSB 
wavelets  have  been  removed  in  these  deconvolutions. 


ANNUAL  REPORT 


Multichannel  Deconvolution  of  P  Waves 


TGAL-86-6 


l 


Results  for  Yucca  Flats  events 

Data  for  Yucca  Rats  events  were  available  at  three  of  the  AWRE  arrays:  EKA,  GBA, 
and  YKA.  The  locations  of  the  events  are  shown  in  Figure  37,  and  deconvolved  source  time 
functions  with  the  VSB  source  pulse  removed  are  shown  in  Figures  38  to  40.  The  recordings 
of  Yucca  Flats  events  exhibit  complex  waveforms  and  ringing  at  all  arrays.  This  is  commonly 
attributed  to  reverberations  in  the  low  velocity  sediments  at  Yucca  Flats  (for  finite  difference 
simulations,  see  McLaughlin  et  al  1986a).  It  appears  that  merely  factoring  out  the  instrument 
and  t*  operators  is  not  sufficient  and  some  deconvolution  of  the  near  source  ringing  would 
have  to  be  introduced.  At  this  point,  this  would  have  gone  beyond  the  immediate  objectives  of 
this  project,  although  we  plan  to  attempt  this  in  the  future.  The  waveforms  at  YKA  are  also 
triplicated  by  upper  mantle  discontinuities.  The  late,  ringing  wavetrain  about  12  seconds  after 
the  first  arrival  in  many  of  the  waveforms  in  Figure  40  is  probably  due  to  the  forward  branch 
associated  with  the  600  km  discontinuity.  Tentative  picks  for  some  "pP"  arrivals  for  Yucca 
Flats  show  them  to  be  quite  late,  but  the  ringing  makes  such  picks  doubtful.  Since  NTS  is  a 
low  Q  site  (Der  et  al  1982,  1985)  we  have  a  very  limited  bandwidth  to  work  with  at 
telesesimic  distances.  The  trace  reconstructions  for  the  Yucca  Rats  events  are  good  reproduc¬ 
tions  of  the  original  traces  (Figure  41).  Often  the  reconstructions  seem  to  have  better  signal- 
to-noise  ratios  than  the  original  traces.  This  is  not  surprising  since  our  deconvolution  method 


can  be  considered  as  an  array  processor  that  uses  the  redundancy  of  the  data  set  to  enhance  the 
source  and  site  factors.  Moreover,  by  limiting  the  frequency  ranges  in  the  deconvolutions  and 
reconstructions  we  have  reduced  the  amplitudes  of  the  low  and  high  frequency  background 
noise  seen  in  the  original  records. 

The  data  set  of  events  recorded  at  EKA  (Figure  38)  is  particularly  interesting  in  that  it 
includes  Piledriver,  one  of  the  few  granite  shots  at  NTS  and  an  important  calibration  event  in 
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Figure  37.  Locatir  >  of  the  Yucca  Flats  events  analyzed  in  this  report. 
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Figure  38.  Source  time  function  estimates  for  a  set  of  Yucca  Flats  events  at  EKA.  VSB 
wavelets  have  been  removed  in  these  deconvolutions.  The  explosion  Piledriver  (north  of 
Yucca  Rats)  she  ws  a  simpler  waveform  without  the  ringing  associated  with  Yucca  Flats 
events. 
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Figure  39.  Source  time  function  estimates  for  a  set  of  Yucca  Flats  events  at  GBA.  VSB 
wavelets  have  been  removed  in  these  deconvolutions. 
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Figure  40.  Source  time  function  estimates  for  Yucca  Flats  events  at  YKA.  The  second  set  ot 
arrivals  is  an  upper  mantle  triplication.  VSB  wavelets  have  been  removed  in  these  deconvolu¬ 
tions. 


Figure  41.  Examples  of  trace  reconstructions  (for  event  Cabrillo)  at  YKA  as  compared  to  ori¬ 
ginal  traces. 
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many  studies.  This  event  is  not  located  at  Yucca  Flats  but  in  a  gianite  body  north  of  it  Due 
to  the  shon  distance  between  Yucca  Flats  and  the  location  of  Piledriver,  the  recording  site  fac¬ 
tors  for  Yucca  Flats  can  also  be  expected  to  be  applicable  for  Piledriver.  The  deconvolution 
for  Piledriver  it  is  relatively  simple  compared  to  those  for  the  Yucca  Flats  events.  Pilednvet 
exhibits  a  double  negative  peak  following  the  first  positive  arrival  (P).  The  first  of  these  is 
timed  appropriately  for  the  expected  pP,  though  we  do  not  have  an  explanation  for  the  second 
negative  peak. 

Results  for  three-component  data 

The  particle  motions  of  teleseismic  body  waves  are  only  in  the  first  approximation  linear 
for  P  and  S.  Even  for  simple  layered  structures  of  the  crust  the  motion  becomes  elliptical  in 
the  most  complex  manner  (Su  and  Dorman  1965)  with  prograde  and  retrograde  motions  in  the 
various  frequency  bands.  For  realistic,  three  dimensional  models  of  the  near-station  geology 
the  motion  becomes  truly  complex  and  three  dimensional  even  for  supposedly  simple,  linear  P 
waves.  While  the  first  few  cycles  of  the  P  wave  are  reasonably  linear-elliptical  in  the  plane  of 
the  azimuth  of  the  arrival,  the  transverse  component  quickly  builds  up  with  time,  approaching 
to  a  respectable  fraction  of  the  vertical  and  radial  components  in  energy.  It  is  a  logical  exten¬ 
sion  of  the  deconvolution  work  on  multiple  sites  to  approach  the  problem  by  treating  each 
component  of  motion  as  a  separate  "site"  in  the  deconvolution  process.  Each  component  will 
have  a  "site  factor"  in  the  frequency  domain  which  characterizes  the  transformation  of  the 
motion  due  to  local  geology,  these  factors  will  depend,  of  course,  on  the  slowness  vector 
(azimuth-distance)  of  the  P  arrivals  and  may  thus  be  used  to  recognize  or  characterize  source 
regions  of  P  waves. 

The  approach  developed  here  is  fundamentally  different  from  practically  all  of  the  3- 
component  processing  work  being  done  recently  (Christofferson  et  al  1986,  Harris  1981) 
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instead  of  assuming  a  simply  polarized  wave  appropriate  to  a  halfspace,  we  plan  to  utilize  the 
particle  motion  as  it  is,  with  all  its  complexities.  A  knowledge  base  of  such  characteristics 
may  be  accumulated  observationally  for  a  large  number  of  source  areas,  without  the  need  for 
understanding  the  detailed  mechanism  of  the  generation  of  complex  particle  motion.  Although 
such  an  approach  may  be  less  practical  presently,  advances  of  computer  technology  will  make 
it  quite  practical  to  store  such  information.  We  believe  that  source  regions  may  be  identified 
quickly  and  automatically  by  recognizing  patterns  of  3-component  site  factors  rather  than  com¬ 
puting  back-azimuths  based  on  simplistic  assumptions  about  polarization  patterns.  Any  sys¬ 
tematic  deviations  between  the  true  azimuths  and  the  back  azimuths  can  be  automatically 
accounted  for  in  such  schemes,  as  long  as  the  master  events  used  to  derive  the  3-component 
"site"  factors  are  accurately  located  (possibly  independently). 

To  test  such  concepts  we  have  assembled  a  data  base  of  four  master  events  at  the  Kazakh 
test  site  as  recorded  at  the  RSTN  station  network.  The  deconvolved  source  time  functions 
appear  to  be  quite  simple  for  these,  mainly  consisting  of  P  and  (presumably)  pP  arrivals  as 
marked  (Figure  42).  The  site  impulse  responses  tend  to  be  the  most  complex  on  the  transverse 
components  building  up  gradually  with  time,  as  expected,  but  simpler  and  impulselike  on  the 
vertical  and  radial  components  (Figure  43).  Nevertheless,  all  the  component  "site"  time  func¬ 
tions  seem  to  be  much  more  complicated  at  RSCP  and  RSSD,  two  sites  the  geological  struc¬ 
ture  under  which  have  been  shown  to  be  complex  by  other  studies  (Owens  et  al  1986).  Some 
raw  data  traces  of  these  events  are  shown  in  Figure  44,  and  opposing  them,  we  show  the  trace 
reconstructions  for  easy  comparison.  The  good  similarity  of  the  original  and  reconstructed 
traces  show  that  the  "site"  approach  works  well  for  three  component  data  despite  the  larger 
variation  of  the  back  azimuths  from  the  source  region  compared  to  seismic  arrays,  and  that  the 
complexities  in  particle  motion  can  be  reliably  characterized  by  "component-site"  transfer  func¬ 
tions. 
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Figure  42.  Source  time  function  estimates  for  four  E-EKTS  events  obtame 
data  set  from  RSTN  three  component  recordings.  VSB  wavelets  have  been  i 
deconvolutions.  All  of  these  show  the  prominent  "pP"  arrival. 
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Figure  43.  Site  time  function  estimates  for  all  components  at  the  five  RSTN  stations.  The 
RSCP  and  RSSD  site  functions  are  considerably  more  complex  than  those  for  RSON,  RSNY 
and  RSNT.  This  may  be  due  to  more  complex  geological  structures  at  the  former  stations. 
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Figure  44.  Original  and  reconstructed  traces  at  three  component  RSTN  stations  for  a  E-EKTS 
event.  These  reconstructions  are  quite  good  considering  the  fact  that  the  RSTN  stations  cover 
a  wide  area. 
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USE  OF  DECONVOLUTION  RESULTS  FOR  YIELD  ESTIMATION 

From  the  above  discussions  it  is  apparent  the  the  deconvolution  method  uses  repeated 
averaging  over  sites  and  events  in  an  iterative  manner,  while  deconvolving  the  source  and  site 
factors  as  well  as  the  factors  (instrument,  t*  and  source  wavelet)  assumed  to  be  known.  The 
basic  idea  is  to  average  waveforms  such  that  all  traces  are  given  an  equal  weight.  In  order  to 
ensure  this  the  whole  process  is  started  by  equalizing  all  trace  amplitudes  by  initial  amplitude 
factors.  This  ensures  that  the  largest  event  does  not  dominate  the  site  factor  estimates  and  that 
the  site  with  the  largest  amplitude  does  not  contaminate  the  source  factor  estimates.  This  way 
the  traces  with  both  the  smallest  and  the  largest  absolute  amplitudes  can  be  reconstructed  with 
about  equal  fidelity  to  the  original  waveforms.  An  inevitable  consequence  of  this  approach  is 
that  the  absolute  signal  amplitudes  become  intractable  through  the  process.  The  idea  of 
waveform  averaging,  an  essential  ingredient  of  this  deconvolution  method,  is  in  basic  conflict 
with  the  preservation  of  the  absolute  signal  amplitudes.  To  preserve  the  approximate  equal 
weighting  of  the  trace  waveforms  inside  the  program  the  t*,  instrument  and  source  wavelet 
amplitudes  are  all  normalized,  rather  than  scaled.  In  order  to  go  back  to  absolute  amplitudes 
we  must  relate  the  signal  amplitudes  or  to  the  deconvolution  results  through  the  absolute 
amplitudes  measured  for  the  set  of  traces  originally  entering  the  procedure. 

A  better  measure  of  the  event  size  may  be  constructed  the  following  way;  by  dividing  out 
the  source  and  site  time  functions  and  instrument  response  in  which  the  first  positive  peak  is 
normalized  to  unity  and  a  t*  operator  that  is  not  normalized  the  result  (properly  limited  in 
bandwidth  to  match  the  band  used  for  the  deconvolution)  will  be  a  simple  impulsive  waveform 
the  log(amplitude)  of  which  may  be  used  as  a  new  "mb"  parameter.  The  deconvolution  of  all 
these  factors  constitutes  a  set  of  corrections  for  site  effects,  pP,  instrument  and  t*  with  the 
absolute  amplitude  information  in  the  original  wave  amplitudes  restored.  The  absolute 
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amplification  factors  due  to  sites  will  still  have  to  be  averaged  out  in  the  absence  of  any  abso¬ 
lute  site  amplitude  reference  levels  in  the  data.  During  the  continuation  of  this  work  we  shall 
attempt  to  test  this  new  magnitude  measure. 
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CONCLUSIONS 

Most  of  the  ambiguities  associated  with  the  pP  estimation  methods  proposed  in  the  past 
are  associated  with  the  fact  that,  for  most  explosions,  the  simple  P+pP  model  is  not  valid. 
This  can  be  verified  by  analyses  of  interevent  spectral  ratios  and  the  results  of  our  multichannel 
deconvolutions.  In  addition,  for  single  channels  the  site  effects  will  also  disguise  any  pP 
arrivals  even  if  they  are  present.  Therefore,  most  results  concerning  pP  using  the  P+pP  model, 
narrow  band  synthetics,  cross-equalization  and  spectral  fitting,  can  be  shown  to  be  invalid  by 
broad-band  analyses. 

*  C-EKTS  explosion  deconvolutions  commonly  show  no  or  a  poorly  defined  pP  arrival. 
The  radiation  from  the  C-EKTS  site  seems  to  have  strong  azimuthal  asymmetries  as  evi¬ 
denced  by  the  poor  reconstructions  of  the  original  traces  following  joint  deconvolutions 
of  data  from  all  arrays. 

*  E-EKTS  explosions  deconvolved  appear  to  have  relatively  simple  source  functions  with  a 
direct  P  followed  by  a  pP  arrival  with  a  time  delay  in  the  0.3-0.6  sec  range. 

*  Novaya  Zemlya  explosions  show  a  complex  late  arrival  instead  of  a  well  recognizable  pP. 
YKA  waveforms  and  deconvolutions  are  complex  compared  to  deconvolutions  at  EKA 
and  GBA. 

*  For  Yucca  Flats  events  we  have  obtained  poor  deconvolution  results  because  of  the 
extremely  limited  bandwidth  of  these  signals  and  the  extreme  ringing  at  that  test  site. 
Piledriver  deconvolutions  at  EKA  show  a  simple  waveform  indicating  two  arrivals  with 
negative  polarity  following  P. 
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*  The  French  nuclear  explosions  in  the  Sahara  gave  complex  source  functions  with  strong 
azimuthal  variations. 

*  Application  of  the  multichannel  deconvolution  technique  to  three  component  recordings 
indicates  that  the  site  response  concept  is  applicable  to  3-dimensional  particle  motion. 
This  indicates  that  this  approach  can  have  application  to  three  component  processing. 

*  Reconstructions  from  the  site  and  source  time  series  and  the  source  wavelets  show  excel¬ 
lent  agreement  with  the  original  traces  including  fine  details  in  the  coda.  The  coda  is  par¬ 
tially  due  to  both  near  site  and  source  scattering. 
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FUTURE  WORK 

During  the  next  year’s  work  on  this  contract  we  shall  continue  our  data  analyses  on  addi¬ 
tional  data  sets;  the  Swedish  network,  NORSAR  will  be  also  utilized  in  addition  to  the  AWRE 
arrays.  More  data  from  NTS  will  be  analyzed  in  addition  to  the  Yucca  Flats  explosions.  For 
Novaya  Zemlya  we  shall  seek  data  representing  both  the  northern  and  southern  tests  areas. 

To  obtain  a  better  understanding  of  the  site  functions  we  shall  analyze  20  earthquakes 
and  other  events  with  azimuths  centered  at  but  different  from  those  of  the  Kazakh  and  NTS 
test  sites.  This  will  enable  us  to  test  the  sensitivity  of  site  transfer  functions  to  deviations  in 
azimuth-slowness  vector).  The  second  year’s  work  will  also  include  the  data  analyses  for 
Sinkiang  explosions. 

By  studying  the  azimuthal  asymmetries  in  the  P  wave  radiation  from  various  test  sites 
and  correlating  them  with  previous  studies  of  strain  release  light  may  be  shed  on  the  possible 
effects  of  such  mechanisms  on  the  short  period  P  wave  radiation  and  the  mb  measurement. 
NTS  and  Kazakh  are  the  most  studied  test  sites  with  respect  to  strain  release  with  numerous 
source  mechanism  solutions  and  F  factor  results  which  will  be  utilized. 

By  deconvolving  the  normalized  reconstructions  from  the  original  traces  as  described 
above  we  shall  test  the  new  energy  measure  for  explosion  generated  P  waves  provided  by  an 
integral  of  the  resulting  simple  pulse.  This  new  magnitude  measure  will  be  compared  to  actual 
yields  (if  known)  or  standard  mb.  The  result  thus  produced  will  be  essentially  corrected  for 
site,  attenuation  and  pP  effects. 
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Department  of  Energy 

ATTN:  Mr.  Max  A.  Koontz  (DP-52) 

International  Security  Affairs 
1000  Independence  Avenue 
Washington,  D.C.  20545 

Lawrence  Livermore  National  Laboratory 
ATTN:  Dr.  J.  Hannon  and  Dr.  M.  Nordyke 
University  of  California 
P.O.  Box  808 
Livermore,  CA  94550 

Los  Alamos  Scientific  Laboratory 
ATTN:  Dr.  K.  Olsen  and  Dr.  T.  Weaver 
P.O.  Box  1663 
Los  Alamos,  NM  87544 

Sandia  Laboratories 
ATTN:  Mr.  P.  Stokes 
Geosciences  Department  1255 
Albuquerque,  NM  87185 

OTHER  GOVERNMENT  AGENCIES 

Central  Intelligence  Agency 
ATTN:  Dr.  L.  Turnbull 
OSI/NED,  Room  5G48 
Washington,  D.C.  20505 


U.S.  Arms  Control  and  Disarmament  Agency 
ATTN:  Dr.  M.  Eimer 

Verification  and  Intelligence  Bureau,  Rm  4953 
Washington,  D.C.  20451 

U.S.  Arms  Control  and  Disarmament  Agency 
ATTN:  Mrs.  M.  Hoinkes 
Multilateral  Affairs  Bureau,  Rm  5499 
Washington,  D.C.  20451 

U.S.  Geological  Survey 

ATTN:  Dr.  T.  Hanks 

National  Earthquake  Research  Center 

345  Middlefield  Road 

Menlo  Park,  CA  94025 

U.S.  Geological  Survey 
ATTN:  Dr.  R.  Masse 
Global  Seismology  Branch 
Box  25046,  Stop  967 
Denver  Federal  Center 
Denver,  CO  80225 

UNIVERSITIES 

Boston  College 
ATTN:  Dr.  A.  Kafka 
Western  Observatory 
381  Concord  Road 
Weston,  MA  02193 

California  Institute  of  Technology 
ATTN:  Dr.  D.  Harkrider 
Seismological  Laboratory 
Pasadena,  C  A  91125 

Columbia  University 
ATTN:  Dr.  L.  Sykes 
Lamont-Doherty  Geological  Observatory 
Palisades,  NY  10964 


Cornell  University 

ATTN:  Dr.  M.  Barazangi 

INSTOC 

Snee  Hall 

Ithaca,  NY  14853 
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Harvard  University 
ATTN:  Dr.  J.  Woodhouse 
Hoffman  Laboratory 
20  Oxford  Street 
Cambridge,  MA  02138 

Massachusetts  Institute  of  Technology  3 

ATTN:  Dr.  S.  Soloman,  Dr.  N.  Toksoz,  and  Dr.  T.  Jordan 
Department  of  Earth  and  Planetary  Sciences 
Cambridge,  MA  02139 

Southern  Methodist  University  1 

ATTN:  Dr.  E.  Herrin 
Geophysical  Laboratory 
Dallas,  TX  75275 

State  University  of  New  York  at  Binghamton  1 

ATTN:  Dr.  F.  Wu 

Department  of  Geological  Sciences 

Vestal,  NY  13901 

Sl  Louis  University  2 

ATTN:  Dr.  O.  Nuttli  and  Dr.  R.  Herrmann 

Department  of  Earth  and  Atmospheric  Sciences 

3507  Laclede 

St.  Louis,  MO  63156 

The  Pennsylvania  State  University  1 

ATTN:  Dr.  S.  Alexander 

Geosciences  Department 

403  Deike  Building 

University  Park,  PA  16802 

University  of  Arizona  1 

ATTN:  Dr.  T.  Wallace 
Department  of  Geosciences 
Tucson,  AZ  85721 

University  of  California,  Berkeley  1 

ATTN:  Dr.  T.  McEvilly 

Department  of  Geology  and  Geophysics 

Berkeley,  CA  94720 

University  of  California  Los  Angeles  1 

ATTN:  Dr.  L.  Knopoff 
405  Hilgard  Avenue 
Los  Angeles,  CA  90024 
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University  of  California,  San  Diego 
ATTN:  Dr.  J.  Orcutt 
Scripps  Institute  of  Oceanography 
La  Jolla,  CA  92093 

University  of  Colorado 
ATTN:  Dr.  C.  Archambeau 
CIRES 

Boulder,  CO  80309 

University  of  Illinois 
ATTN:  Dr.  S.  Grand 
Department  of  Geology 
1301  West  Green  Street 
Urbana,  IL  61801 

University  of  Michigan 
ATTN:  Dr.  T.  Lay 
Department  of  Geological  Sciences 
Ann  Arbor,  MI  48109-1063 

University  of  Nevada 
ATTN:  Dr.  K.  Priesdey 
Mackay  School  of  Mines 
Reno,  NV  89557 

University  of  Southern  California 
ATTN:  Dr.  K.  Aki 
Center  for  Earth  Sciences 
University  Park 
Los  Angeles,  CA  90089-0741 

DEPARTMENT  OF  DEFENSE  CONTRACTORS 

Applied  Theory,  Inc. 

ATTN:  Dr.  J.  Trulio 
930  South  La  Brea  Avenue 
Suite  2 

Los  Angeles,  CA  90036 

Center  for  Seismic  Studies 
ATTN:  Dr.  C.  Romney  and  Mr.  R.  Perez 
1300  N.  17th  Street,  Suite  1450 
Arlington,  VA  22209 

ENSCO,  Inc. 

ATTN:  Mr.  G.  Young 
5400  Port  Royal  Road 
Springfield,  VA  22151 


ENSCO,  Inc. 

ATTN:  Dr.  R.  Kemerait 
445  Pineda  Court 
Melbourne,  FL  32940 

Gould  Inc. 

ATTN:  Mr.  R.  J.  Woodard 
Chesapeake  Instrument  Division 
6711  Baymeado  Drive 
Glen  Burnie,  MD  21061 

Pacific  Sierra  Research  Coip. 

ATTN:  Mr.  F.  Thomas 
12340  Santa  Monica  Boulevard 
Los  Angeles,  CA  90025 

Rockwell  International 
ATTN:  B.  Tittmann 
1049  Camino  Dos  Rios 
Thousand  Oaks,  CA  91360 

Rondout  Associates,  Inc. 

ATTN:  Dr.  P.  Pomeroy 

P.O.  Box  224 

Stone  Ridge,  NY  12484 

Science  Applications,  Inc. 

ATTN:  Dr.  T.  Bache,  Jr. 

P.O.  Box  2351 
La  Jolla,  CA  92038 

Science  Horizons 

ATTN:  Dr.  T.  Cherry  and  Dr.  J.  Minster 
710  Encinitas  Blvd. 

Suite  101 

Encinitas,  CA  92024 

Sierra  Geophysics,  Inc. 

ATTN:  Dr.  R.  Hart  and  Dr.  G.  Mellman 
11255  Kirkland  Way 
Kirkland,  WA  98124 

SRI  International 
ATTN:  Dr.  A.  Florence 
333  Ravensworth  Avenue 
Menlo  Park,  CA  94025 

S-Cubed,  A  Division  of 
Maxwell  Laboratories  Inc. 

ATTN:  Dr.  S.  Day 
P.O.  Box  1620 
La  Jolla,  CA  92038-1620 
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S-Cubed,  A  Division  of 
Maxwell  Laboratories  Inc. 

ATTN:  Mr.  J.  Murphy 
11800  Sunrise  Valley  Drive 
Suite  1212 
Reston,  VA  22091 

Teledyne  Geotech 

ATTN:  Dr.  Z.  Der  and  Mr.  W.  Rivers 
314  Montgomery  Street 
Alexandria,  VA  22314 

Woodward-Clyde  Consultants 

ATTN:  Dr.  L.  Burdick  and  Dr.  J.  Barker 

556  El  Dorado  St. 

Pasadena,  C A  91105 


NON-U. S.  RECIPIENTS 


National  Defense  Research  Institute 
ATTN:  Dr.  O.  Dahlman 
Stockholm  80,  Sweden 

Blacknest  Seismological  Center 

ATTN:  Mr.  P.  Marshall 

Atomic  Weapons  Research  Establishment 

UK  Ministry  of  Defence 

Brimpton,  Reading  RG7-4RS 

United  Kingdom 

NTNF  NORSAR 
ATTN:  Dr.  F.  Ringdal 
P.O.  Box  51 
N-2007  Kjeller 
Norway 


OTHER  DISTRIBUTION 


To  be  determined  by  the  project  office 


